title: GIS Assignment 2: External Data Sources
packages<-c("cowplot","dismo","leaflet","maps","mapdata","raster","rasterVis","readxl","rgbif","rgdal","tidyverse","utils","rdryad")
sapply(packages, library, character.only=T)
library(ggplot2)
library(maps)
library(OpenStreetMap)
library(mapproj)
library(ggmap)
## Google's Terms of Service: https://cloud.google.com/maps-platform/terms/.
## Please cite ggmap if you use it! See citation("ggmap") for details.
## 
## Attaching package: 'ggmap'
## The following object is masked from 'package:dismo':
## 
##     geocode
## The following object is masked from 'package:cowplot':
## 
##     theme_nothing
chrom_dismo <- gbif("chrosomus", species = "erythrogaster", ext = c(-91,-81, 36, 40),
                   geo = TRUE, sp = TRUE, download = TRUE,
                   removeZeros = TRUE)
chrom_dismo_df <- cbind.data.frame(chrom_dismo@coords[,1],
                                  chrom_dismo@coords[,2])
colnames(chrom_dismo_df) <- c("x","y")
us <- map_data("state")
ggplot(data = chrom_dismo_df, aes(x=x, y=y)) +
  geom_polygon(data = us, aes(x=long, y = lat, group = group),
               fill = "white", color="black") +
  geom_point() + xlab("Longitude") + ylab("Latitude") +
  coord_fixed(xlim = c(-93,-76), ylim = c(33, 60)) +
  xlab("Longitude") + ylab("Latitude") + ggtitle("Chromsomus erythrogaster") + 
  theme_bw() + theme(plot.title = element_text(hjust = 0.5)) + 
  theme(panel.grid.major = element_blank(), 
        panel.grid.minor = element_blank(),
        panel.background = element_rect(fill = "lightblue"))

Here I limited my species counts to 2000 individuals per species

chrom_rgbif <- occ_data(scientificName = "Chrosomus erythrogaster",
                 hasCoordinate = TRUE, limit = 2000,
                 decimalLongitude = "-130, -40", 
                 decimalLatitude = "30, 80")

cumb_rgbif <- occ_data(scientificName = "Chrosomus cumberlandensis",
                 hasCoordinate = TRUE, limit = 2000,
                 decimalLongitude = "-130, -40", 
                 decimalLatitude = "30, 80")
Figure 2. Blackside Dace (Chrosomus cumberlandensis)
Figure 3. Southern Redbelly Dace (Chrosomus erythrogaster)
chrom_rgbif_df <- cbind.data.frame(chrom_rgbif$data$species,
                                  chrom_rgbif$data$decimalLatitude,
                                  chrom_rgbif$data$decimalLongitude,
                                  chrom_rgbif$data$stateProvince,
                                  chrom_rgbif$data$year)


cumb_rgbif_df <- cbind.data.frame(cumb_rgbif$data$species,
                                  cumb_rgbif$data$decimalLatitude,
                                  cumb_rgbif$data$decimalLongitude,
                                  cumb_rgbif$data$stateProvince,
                                  cumb_rgbif$data$year)
colnames(chrom_rgbif_df) <- c("species","y","x","state","year")
colnames(cumb_rgbif_df) <- c("species","y","x","state","year")
chrom_rgbif_df <- chrom_rgbif_df[complete.cases(chrom_rgbif_df[1:4]),]
cumb_rgbif_df <- cumb_rgbif_df[complete.cases(cumb_rgbif_df[1:4]),]
ggplot() +
  geom_polygon(data = us, aes(x=long, y = lat, group = group),
               fill = "white", color="black") +
  geom_point(data = chrom_rgbif_df, aes(x=x, y=y, color = species), size = 3) +
  geom_point(data = cumb_rgbif_df, aes(x=x, y=y, color = species), size = 3) +  
  coord_fixed(xlim = c(-107,-65), ylim = c(24,50)) +
  xlab("Longitude") + ylab("Latitude") + ggtitle("Chrosomus erythrogaster and Chrosomus cumberlandensis") + 
  guides(color=guide_legend("Legend", override.aes = list(size = 4))) +
  theme_bw() + theme(plot.title = element_text(hjust = 0.5)) + 
  theme(legend.position = "bottom") +
  theme(legend.title.align = 0.5, legend.box.just = "center") +
  theme(panel.grid.major = element_blank(), 
        panel.grid.minor = element_blank(),
        panel.background = element_rect(fill = "lightblue"))

pal <- colorNumeric(
  palette = "RdYlBu",
  domain = chrom_rgbif_df$year)
colors <- colorFactor(c("#F8766D","#00BA38","#619CFF"), chrom_rgbif_df$year)
leaflet(chrom_rgbif_df) %>% 
  addTiles() %>% 
  addCircleMarkers(chrom_rgbif_df$x,
                   chrom_rgbif_df$y,
                   popup = chrom_rgbif_df$species,
                   label= chrom_rgbif_df$year,
                   color = ~pal(year),
                   weight = 1,
                   fillOpacity = 0.7) %>%
  addMiniMap(position = 'topright',
             width = 100, 
             height = 100,
             toggleDisplay = FALSE) %>%
  addScaleBar(position = "bottomright")%>%
  addLegend("bottomleft", pal = pal, values = chrom_rgbif_df$year,
    title = "Year Recorded",
    labFormat = labelFormat(),
    opacity = 1)
Figure 1. Dung Fly
Dung_Fly <- read.csv("Dung_Fly.csv")
View(Dung_Fly)
canada <- getData('GADM' , country = "CAN", level=0)
world <- map_data("worldHires")
main_map <- ggplot(Dung_Fly, aes(Longitude, Latitude)) +
  geom_polygon(data = us, aes(x=long, y = lat, group = group),
               fill = "gray", color="white") +
   geom_polygon(data = canada, aes(x=long, y = lat, group = group),
               fill = "gray", color="white") +
  geom_point(aes(color = Altitude, size = `Clutch.Size..First.clutch.`)) +
  coord_fixed(xlim = c(-125,-66), ylim = c(20,60)) +
  xlab("Longitude") + ylab("Latitude") + ggtitle("Dung Fly Clutch Size and Altitude") + 
  guides(color=guide_legend("Altitude", override.aes = list(size = 3))) + 
  guides(size=guide_legend("Clutch.Size..First.clutch.")) +  
  theme_bw() + theme(plot.title = element_text(hjust = 0.5)) + 
  theme(legend.position = "right") +
  theme(legend.title.align = 0.5, legend.box.just = "center") +
  theme(panel.grid.major = element_blank(), 
        panel.grid.minor = element_blank(),
        panel.background = element_rect(fill = "lightblue"))
main_map

bioclim <- getData(name = "worldclim", res = 2.5, var = "bio", path = "./")
names(bioclim) <- c("Ann Mean Temp","Mean Diurnal Range","Isothermality","Temperature Seasonality",
                    "Max Temp Warmest Mo","Min Temp Coldest Mo","Ann Temp Range","Mean Temp Wettest Qtr",
                    "Mean Temp Driest Qtr","Mean Temp Warmest Qtr","Mean Temp Coldest Qtr","Annual Precip",
                    "Precip Wettest Mo","Precip Driest Mo","Precip Seasonality","Precip Wettest Qtr",
                    "Precip Driest Qtr","Precip Warmest Qtr","Precip Coldest Qtr")
bio_extent <- extent(x = c(
  min(chrom_rgbif_df$x),
  max(chrom_rgbif_df$x),
  min(chrom_rgbif_df$y),
  max(chrom_rgbif_df$y)))
bioclim_extent <- crop(x = bioclim, y = bio_extent)
bioclim_model <- bioclim(x = bioclim_extent, p = cbind(chrom_rgbif_df$x,chrom_rgbif_df$y))
presence_model <- dismo::predict(object = bioclim_model, 
                                 x = bioclim_extent, 
                                 ext = bio_extent)
gplot(presence_model) + 
  geom_polygon(data = us, aes(x= long, y = lat, group = group),
               fill = "gray", color="black") +
  geom_raster(aes(fill=value)) +
  geom_polygon(data = us, aes(x= long, y = lat, group = group),
               fill = NA, color="black") +
  geom_point(data = chrom_rgbif_df, aes(x = x, y = y), size = 2, color = "black", alpha = 0.5) +
  scale_fill_gradientn(colours=c("brown","yellow","darkgreen"), "Probability") +
  coord_fixed(xlim = c(-107,-65), ylim = c(24,50)) +
  xlab("Longitude") + ylab("Latitude") + ggtitle("Probability of PILO Occurrence") + 
  theme_bw() + theme(plot.title = element_text(hjust = 0.5)) + theme(legend.position = "right") +
  theme(panel.grid.major = element_blank(), 
        panel.grid.minor = element_blank(),
        panel.background = element_rect(fill = "lightblue"))

colors <- c("brown","yellow","darkgreen")
leaflet() %>% 
  addTiles() %>%
  addRasterImage(presence_model, colors = colors, opacity = 0.8) %>%
  addCircleMarkers(chrom_rgbif_df$x,
                   chrom_rgbif_df$y,
                   weight = 1,
                   color = "green",
                   fillColor = "blue",
                   fillOpacity = 0.7) %>%
  addMiniMap(position = 'topright',
             width = 100, 
             height = 100,
             toggleDisplay = FALSE) %>%
  addScaleBar(position = "bottomright")
state <- map_data("state")
county <- map_data("county")
point <- data.frame("x" = -84.7471, "y" = 34.6953)
View(chrom_rgbif_df)
inset <- ggplot() + 
  geom_polygon(data = us, aes(x=long, y = lat, group = group),
               fill = "yellow", color="green") +
  geom_polygon(data = canada, aes(x=long, y = lat, group = group),
               fill = "pink", color="green") +
  geom_point(data = Dung_Fly, aes(x=Longitude, y=Latitude, color = Clutch.Size..First.clutch.), size = 5) +
  coord_map(xlim = c(-130, -70), ylim = c(24,50), "polyconic") + 
  theme(panel.background = element_rect(fill = "lightblue"), panel.grid.major = element_blank(), panel.grid.minor = element_blank(), 
        axis.line=element_blank(), axis.text.x=element_blank(), axis.text.y=element_blank(),axis.ticks=element_blank(), 
        axis.title.x=element_blank(), axis.title.y=element_blank()) +
theme(plot.margin=grid::unit(c(0,0,0,0), "mm"))
inset

ggdraw() +
draw_plot(main_map) + 
draw_plot(inset, x = 0.50, y = 0.50, width = 0.50, height = 0.50)